Chytrid infections exhibit historical spread and contemporary seasonality in a declining stream-breeding frog

Species with extensive geographical ranges pose special challenges to assessing drivers of wildlife disease, necessitating collaborative and large-scale analyses. The imperilled foothill yellow-legged frog (Rana boylii) inhabits a wide geographical range and variable conditions in rivers of California and Oregon (USA), and is considered threatened by the pathogen Batrachochytrium dendrobatidis (Bd). To assess drivers of Bd infections over time and space, we compiled over 2000 datapoints from R. boylii museum specimens (collected 1897–2005) and field samples (2005–2021) spanning 9° of latitude. We observed a south-to-north spread of Bd detections beginning in the 1940s and increase in prevalence from the 1940s to 1970s, coinciding with extirpation from southern latitudes. We detected eight high-prevalence geographical clusters through time that span the species' geographical range. Field-sampled male R. boylii exhibited the highest prevalence, and juveniles sampled in autumn exhibited the highest loads. Bd infection risk was highest in lower elevation rain-dominated watersheds, and with cool temperatures and low stream-flow conditions at the end of the dry season. Through a holistic assessment of relationships between infection risk, geographical context and time, we identify the locations and time periods where Bd mitigation and monitoring will be critical for conservation of this imperilled species.


Introduction
Threatened species with large geographical ranges often require unique, regional conservation strategies to combat stressors such as infectious disease.Pathogen surveys and reporting have become standard for North American wildlife diseases [1,2]; however, relative risk across a landscape and among populations within species remains difficult to anticipate, especially when data are collected by separate research groups [3].Central reporting databases [4], synthetic analyses and retrospective surveys can help assess disease threats and identify high-risk populations.
Among the most significant wildlife diseases, amphibian chytridiomycosis caused by the fungal pathogen Batrachochytrium dendrobatidis (Bd) has contributed to declines of hundreds of species worldwide [5]; but see [6].In North America, notable Bd-associated declines have occurred across the west including the southern Rocky Mountains [7,8], Arizona and New Mexico [9,10], Nevada [11] and California [12][13][14].In several of these cases, infection outcomes varied widely among populations due to host-related and environmental factors including genetics, prior Bd exposure and abiotic conditions [15][16][17].
For the stream-dwelling foothill yellow-legged frog, Rana boylii, Bd's role in the species' changing abundance across its endemic range (California and Oregon, USA) is not well-understood.The species has declined for at least the last half-century, with extirpations reported from xeric lower latitudes [18], at the wetter northern range limit [19] and downstream of large dams range-wide [20].A mix of abiotic and biotic factors influence Bd infection risk and disease dynamics in many systems, including elevation, latitude, climate, habitat quality and host characteristics [21].The relative importance of these factors remains unclear in rivers with winter flood/summer drought flow regimes typical across R. boylii's geographical range.Bd is considered a significant potential threat to R. boylii [22] because it is implicated in the species' disappearance from rivers of California's South Coast [23] and in recent autumn die-offs of R. boylii in Central Coast streams [24,25].A large-scale assessment of Bd infections is needed to clarify how infections relate to historical declines in some regions' rivers and persistence in others, identify clusters of increased infection risk across the species' range, and evaluate how infection incidence and severity changes with the seasonality of the Mediterranean climate and across the diverse ecoregions that R. boylii occupies.
Here, we leverage data from over 2000 field and museum samples covering 124 years to synthesize knowledge and evaluate patterns of Bd infections in R. boylii.We use a combination of modelling approaches and spatial scan statistics to ask: (i) how are Bd detections in R. boylii are distributed over space and time, (ii) whether watersheds with high versus low Bd infection risk clustered historically and today, and (iii) how Bd infections are related to biotic and abiotic factors.Our results highlight priority populations for Bd mitigation, regions that are data-deficient and warrant further sampling and monitoring, and remaining gaps in our knowledge about Bd susceptibility in R. boylii.Our study serves as a resource for wildlife managers implementing disease mitigation and species recovery projects, such as re-introductions, and as an example of collaborative research to address conservation challenges in wide-ranging imperilled species.

Study species
Rana boylii is a stream-dwelling frog (38-81 mm snout-vent length, SVL) [26] whose life cycle timing and breeding migrations between small tributaries and larger main-stem channels occur in synchrony with the winter flood/summer drought flow regime typical of Mediterranean climate rivers.During the transition from high to low flows, adults congregate at lek sites where males reside for weeks to call from underwater, and females visit long enough to mate and attach eggs to rocks in flow-protected locations.Tadpoles reach metamorphosis three-four months later, well before the return of high-flow disturbance the next winter.
These frogs were once common in the rivers and streams flowing through a variety of biomes including moist coniferous forests, oak-savannahs, chaparral and deserts in California and Oregon.Historically, R. boylii ranged across at least 12°of latitude, from the Willamette drainage in Oregon to at least the San Gabriel drainage in southern California [26,27], and from sea level to approximately 1500 m [28].The species tolerates a wide range of conditions, all of which can influence population size and dynamics, including: precipitation that ranges south-to-north from less than 50 to greater than 300 cm yr −1 [29]; seasonally variable stream-flow [30]; primary productivity gradients [31]; and thermal regime gradients from rain-dominated coastal watersheds to inland montane (Sierran) snowmelt-driven watersheds [32].Rana boylii comprises six deeply divergent clades that vary in demographic histories, genetic diversity, and California and federal Endangered Species Act protections [33][34][35][36] (figure 1a).royalsocietypublishing.org/journal/rsos R. Soc.Open Sci.11: 231270

Data compilation and visualization
To test for the effects of space, time and abiotic and biotic factors on Bd infections in R. boylii, we gathered georeferenced Bd data from skin swabs and/or tissues collected from museum specimens or fieldsampled frogs from 22 sources (electronic supplementary material, tables S1 and S2).Skin swabs were collected and DNA was extracted from swabs according to a standard protocol [37] with the exception of one museum study that used a DNA extraction method for formalin-fixed specimens [23].Extracts were analysed with either endpoint or quantitative PCR (qPCR) to detect Bd DNA [38].Two fieldcollected samples [39] were analysed using histological examination of skin at the USGS National Wildlife Health Center (Madison, WI).Museum samples were assessed using histology [40], qPCR [12], or a combination of these two techniques [23,24].Bd assay results and detection method, frog life stage, sex and size were included when available (see electronic supplementary material, methods for additional details on data gathering and quality control).The final dataset included 2145 samples, comprising 461 museum samples of post-metamorphic frogs collected 1897-2005 and 1684 field samples (larval and post-metamorphic frogs) collected 2005-2021.A subset of field-sampled frogs (n = 86) were passive integrated transponder (PIT)-tagged and recaptured; a single sample per year was randomly selected for each resampled frog for inclusion in analyses resulting in 2074 total samples and 1613 field samples after filtering (electronic supplementary material, methods).Sampling areas included streams and rivers in 173 sub-watersheds (12-digit USGS hydrologic unit code, HUC-12) across California and Oregon; the diversity of the watersheds is manifest in their spanning the ancestral homelands of at least 41 Indigenous Peoples past, present and future (figure 1a, electronic supplementary material, tables S3 and S4).

Analyses of Batrachochytrium dendrobatidis detection, prevalence and load over space and time
We tested for a spatio-temporal pattern of Bd detections in museum and field samples using a generalized additive model (GAM; [41]).For the GAM, we grouped samples by HUC-12 and capture year, and filtered the data to only Bd-positive samples.To account for uneven sample sizes across space and time (e.g.no samples collected from southern localities after extirpations in the 1970s), we included total samples (Bd+ and Bd−, n = 2074 after filtering recaptures) at each HUC-12 and capture year combination as a fixed effect.We applied a thin plate spline smoothing function to capture year [42].The GAM was run with the following formula: latitude∼total samples + smoothed capture year.
To identify geographical clusters of Bd infection risk, we used circular spatial scan statistics implemented in SaTScan [43,44].SaTScan identifies geographical clusters of cases (infections) by calculating relative risk (RR) based on the incidence of infections inside versus outside circular windows.Window size is selected using maximum likelihood, and clusters are designated high or low rate based on RR values (RR > 1 indicating more infections than expected in high-rate clusters, RR < 1 indicating fewer infections than expected in low-rate clusters).We performed separate SaTScan analyses for museum versus field samples.For both we applied SaTScan's Bernoulli model to Bd status (Bd±; [44]).We aggregated Bd detection data at the HUC-12 level to improve cluster detection power in areas with fewer samples [45]), and set maximum window radius to 50 km (the approximate mean HUC-12 size in California).We used minimum high-rate cluster case thresholds of n = 3 for museum samples and n = 10 for field samples to balance different limitations on detection in formalinfixed versus field-collected samples with achieving sufficient statistical power in the models (electronic supplementary material, table S5).For cluster selection, we used the following SaTScan settings: Monte Carlo with 999 replications, adjustment for more likely clusters, 10 iterations and alpha = 0.05.We used SaTScan cluster ID and cluster rate (high versus low) in models described below.

Evaluation of Batrachochytrium dendrobatidis's relationships with environmental and hostassociated factors
We used generalized linear mixed models (GLMMs) in glmmTMB [46]) to analyze the effects of abiotic and biotic factors on three response variables: Bd detection (Bd+/Bd-) in (i) museum and (ii) field samples, and (iii) infection load (Bd internal transcribed sequence (ITS) copies) in field samples.We first constructed intercept-only models for each of the three response variables with/without random effects of SaTScan cluster ID and cluster rate (0 = not in a cluster, 1 = low-rate cluster, 2 = high-rate royalsocietypublishing.org/journal/rsos R. Soc.Open Sci.11: 231270 cluster).We then selected the best intercept-only model as a base for model building with additional covariates.
We used public databases to compile 85 potentially relevant environmental covariates at the HUC-12 level including hydroclimatic and geographical factors.To identify redundant covariates, we calculated Pearson's r for every pair of variables; we selected a single variable from pairs with r ≥ 0.7 (electronic supplementary material, tables S6 and S7, and methods).
For the analysis of Bd detection in museum samples, we included only day length, the 24-month standardized precipitation evapotranspiration index (SPEI; a measure of prolonged drought in which negative values indicate more extreme drought and positive values indicate wetter years; hereafter referred to as 'aridity'), latitude, and decade as covariates; we did not include other covariates as these are based on contemporary measurements and do not reflect the time period of museum collections.We explored season and aridity as interactions with decade, as these factors could not have affected Bd status prior to the earliest Bd detections.We also included interactions with latitude to evaluate how drivers of infection vary in strength regionally.
For the analyses of Bd detection and Bd load in field samples, we first selected the most informative set of non-redundant environmental covariates using a boosted regression tree (BRT) analysis using xgboost in R (v. 3.6.2;[41,47,48]).BRT is well-suited for large complex ecological data as it does not assume normality or linearity; it ignores non-informative predictors and is typically unaffected by outliers; and it can accept numeric, categorical, or binary data [49,50].Seven covariates surpassed a BRT threshold of approximately 5% relative importance: T max , latitude, day length, % agriculture on hydric soil, % wetlands remaining, aridity and elevation (electronic supplementary material, figure S1).These covariates were included in the two field sample models.
While spatial and/or temporal variation in environmental covariates are likely, capture year was highly confounded with location in field samples.Therefore, we included interactions between environmental covariates and latitude but not year in the field sample models.We included SVL in the Bd load model based on previously demonstrated relationships [51].We used log 10 -transformed ITS values for Bd load model analyses.
All covariates were centred and scaled for model building.Alternative models were constructed for each response variable with and without interactions.We used Akaike information criterion (AIC) to select the best model for each response variable.

Spatial and temporal distribution of Batrachochytrium dendrobatidis detections
Sampling was geographically and temporally uneven.Museum specimens made up the majority of South Coast and Southern Sierra samples.Field samples from North Coast and Central Coast clades were concentrated in watersheds that were the focus of research projects or have frog population monitoring programmes stipulated in dam operating licenses (electronic supplementary material, table S3 and figure S2).

Abiotic correlates of Batrachochytrium dendrobatidis infections
The best-fit GLMM for museum sample Bd detections included cluster rate as a random intercept; decade as a continuous fixed effect; and one-way interactions between decade and three environmental royalsocietypublishing.org/journal/rsos R. Soc.Open Sci.11: 231270 predictors: day length, aridity and latitude (electronic supplementary material, table S8).Decade was the only significant predictor of Bd detection probability, with a stronger positive relationship between Bd detection probability and time in high-rate clusters relative to samples outside of clusters (decade conditional estimate = 1.54, p < 0.0001; electronic supplementary material, figure S9).The best GLMM for field sample Bd detections included cluster ID and cluster rate as random intercepts; all variables selected by BRT analysis (elevation, day length, latitude, T max , aridity, % agriculture on hydric soils and % wetlands remaining) as fixed effects; and interactions between latitude and two climatic predictors, T max and aridity (electronic supplementary material, table S9).Bd detection probability increased with shorter day length (conditional estimate = −0.16,p < 0.05), lower T max (i.e, cooler temperatures; conditional estimate = −0.62,p < 0.0001) and lower aridity (higher values of SPEI index); conditional estimate = 0.21, p < 0.01).The effects of day length, T max and aridity on Bd detection probability were marginally stronger in high-rate clusters (figure 4a,b, electronic supplementary material, figure S10).Relationships between Bd detection probability and T max varied across latitude (latitude-T max interaction: conditional estimate = 0.36, p < 0.0001).The aridity index showed opposite relationships with observed Bd prevalence in Coastal and Sierran clades: in Coastal clades the relationship was negative (i.e. higher Bd prevalence during prolonged droughts; simple linear regression (SLR), t = −10.05,p < 0.01, R 2 = 0.06), while in Sierran clades it was positive (higher Bd prevalence outside of droughts; SLR, t = 12.2, p < 0.05, R 2 = 0.06; figure 4c,d).Observed Bd prevalence showed a significant negative relationship with T max in Coastal clades only (SLR, t = −4.15,p < 0.001, R 2 = 0.23; figure 4e,f).
The best model for field sample Bd load included cluster ID as a random intercept; T max , day length, aridity, % agriculture on hydric soils, % wetlands remaining, latitude, elevation and SVL as fixed effects; and no interactions between latitude and other predictors (electronic supplementary material, table S10).Higher Bd loads were associated with lower aridity (conditional estimate = 0.51, p < 0.01), smaller SVL (conditional estimate = −0.85,p < 0.0001), and lower T max (conditional estimate = −0.35,p < 0.05).

Spatio-temporal patterns in Batrachochytrium dendrobatidis infections
We analysed spatio-temporal patterns in Bd detections in R. boylii using a compiled dataset comprising over 2000 samples spanning 124 years.Our analyses revealed several notable patterns.First, the 39 Bd+ museum specimens coupled with the earliest Bd detections from the field show a south-to-north progression of detections until 2005 when field swabbing began.Second, Bd prevalence increased through the decades, but the trajectory and timing of peak prevalence varied among climate zones.Among extant populations, those from the most xeric region, California's Central Coast (Central Coast and South Coast clades) had the highest Bd prevalence and infection loads; this region also had some of the earliest Bd detections historically.Third, observed Bd infections varied by frog sex and life stage: Bd prevalence was higher in adult males than females, while loads were highest in juveniles in autumn ( primarily from the Central Coast clade).Lastly, we found a strong influence of Mediterranean climate seasonality, whereby both prevalence and loads were lowest during warm and dry summer months.Below we discuss potential implications of these findings, limitations of our analyses and suggestions for future work.

Batrachochytrium dendrobatidis dispersal across the US west coast
The earliest Bd detections in R. boylii across California and Oregon exhibited a pattern of northward spread corroborated by three historical high-rate clusters: a South Coast cluster in the 1960s-1970s royalsocietypublishing.org/journal/rsos R. Soc.Open Sci.11: 231270 (San Gabriel watershed), a Central Coast cluster in the 1960s-1980s (San Francisco Bay Area watersheds), and a Northern Sierra/Feather River cluster from the 1970s to 2000s (Yuba and Feather watersheds, Sacramento River drainage).Bd prevalence also increased across the decades, rising from less than 5% between 1890 and 1969 to greater than 50% by the 1970s, and staying moderately high (approx.30%) until the 2000s.Moreover, decade was the only significant predictor of Bd prevalence in museum samples.Taken together, our results suggest that Bd infections in R. boylii coincided with historical declines in the South Coast of California before spreading northward across California and Oregon during the middle of the twentieth century.These findings are consistent with previous multi-species museum studies reporting that Bd spread in a wave-like pattern across Central and South America [52][53][54] and northward in southern California [24,55], and that Bd prevalence increased in California in the 1960s-1970s [12,56,57].Though our results may suggest directional spread of Bd in R. boylii, the 'Bd wave hypothesis' may not fully explain the history of chytridiomycosis in our sampling area and neighbouring regions.Previous retrospective studies detected Bd in museum specimens collected decades before reports of chytridiomycosis-associated declines from Baja California, Mexico (1894; [58]), southern California (1915; [23]) and Costa Rica (1964; [57]).Additional evidence suggests that Bd may have been introduced multiple times in California [59] and may have spread from human population centres including the San Francisco Bay Area [40] and the Los Angeles Basin [23].Time-dated phylogenetic analyses also suggest that Bd may have been present in California for approximately 500-1,400 years before present; however, rapid evolution at high elevation sites may influence these estimates [60].

Batrachochytrium dendrobatidis prevalence and load depends on sex and life stage
Behaviour and phenology may explain infection differences among sexes and life stages.Higher prevalence in males may be due to congregating at lek (group breeding) sites and territorial wrestling [61,62] that could increase Bd transmission, as well as reproductive endocrines and stress hormones that can mediate immune function variation among sexes [63][64][65].Substantial time spent calling underwater [66] could also increase Bd growth on males, as the fungus is aquatic and psychrophilic (cold-preferring).Water temperatures when R. boylii breed are typically cool, averaging 10-12°C [67,68].Whether males experience increased negative effects of infections is unknown, but in four other California and Oregon ranids, males exhibit greater Bd-associated mortality than females [69].The observation that juvenile R. boylii had the highest pathogen loads may be attributed to the fact that life stage and season are inextricably coupled.Young of the year are abundant in cooler late summer and autumn months when environmental temperatures are similar to the thermal growth optimum for Bd [31,32,[70][71][72].

A Batrachochytrium dendrobatidis hotspot in California's Central Coast
Two Central Coast watersheds, Alameda and Coyote Creeks, grouped into a high-rate cluster and exhibited the highest Bd loads.Bd sampling at these sites was triggered by observations of sick and dying frogs in autumn; thus these are not random samples.Nonetheless, factors specific to these watersheds may contribute to greater chytridiomycosis mortality, such as intermittent surface flow during dry seasons that crowds frogs at remnant pools and promotes Bd transmission; frequent human visitation facilitating transmission of multiple Bd strains; and the presence of non-native reservoir/amplification hosts (bullfrogs and crayfish; [24,25,73]).Moreover, the Central Coast clade has especially low genetic diversity [35], which may increase chytridiomycosis susceptibility [74,75].

Do season, temperature and hydrologic conditions drive Batrachochytrium dendrobatidis infections?
Shorter day length and lower monthly maximum temperatures were associated with higher Bd infection probability and loads.Our models also indicated Bd infections were associated with successive wet years, but the influence of temperature and aridity varied significantly among clades and latitudes.Moreover, T max and aridity serve as reasonable proxies for water temperature and stream-flow, but have limitations.Air and water temperatures in rivers are highly correlated [76], but can be decoupled by factors including stream-bed geomorphology, canopy cover and hypolimnetic releases from upstream impoundments [77][78][79][80].Stream-flows should decrease during droughts, but can be royalsocietypublishing.org/journal/rsosR. Soc.Open Sci.11: 231270 maintained by dam releases [20], intensified by water abstraction [81] and vary with hillslope water storage capacity [82].Given the significant site-to-site variation in these factors affecting water temperature and stream-flow that could mask the associations of Bd with T max and aridity, we interpret the overall relationships as being valid but requiring consideration of local context.Increasing Bd infections and loads with decreasing temperatures and greater stream-flow is consistent with previous reports of seasonal and latitudinal influences on Bd infection in North America (e.g.[83]) and Bd's low thermal preferences [84].Yet, the effects of temperature and aridity were region-dependent.Snowmelt-driven Sierran watersheds exhibited higher Bd prevalence in 2017-2018 when precipitation was higher, whereas prevalence was higher in rain-fed coastal watersheds during prolonged drought.Coastal samples were collected during more extreme drought conditions (SPEI < −2) compared with Sierran watersheds (min SPEI = −1), when frogs concentrate in shrinking pools during the cool dry autumn, probably increasing Bd transmission [24,25].Extreme drought is expected to intensify with climate change across R. boylii's range [85].

Persisting with a deadly pathogen
The range-wide ubiquity of Bd we document here and episodes of chytridiomycosis related die-offs in Central Coast populations [25] raise questions about whether the persistence of robust R. boylii populations on the North Coast of California and almost complete extirpation from xeric regions of the South Sierra and South Coast are attributed to differences in Bd susceptibility and/or environmental drivers.Mechanisms associated with variation in Bd susceptibility on the individual level include host genetics; innate, acquired or trans-generational ability to clear infections; anti-Bd function of skin secretions and/or microbiome; and microclimate [74,75,[86][87][88].Antimicrobial skin peptides (AMPs) have been characterized in R. boylii [89,90], but the relationships between AMPs and Bd infection remain to be explored in this species.Cases among the set of recaptured frogs in which infection had cleared (n = 18) or load was substantially reduced (n = 9) upon recapture 1-3 years later could be due to any of the mechanisms listed above.
Behaviour and environmental conditions probably also play roles in the ability to clear Bd infections.Rana boylii is known to bask on river banks and emergent rocks in the stream channel [26,91].Short-term elevated temperature treatments that mimic basking can inhibit Bd growth [92,93].Warm and precipitation-free summers with long days provide ample basking opportunities across R. boylii's range, potentially explaining the low summer Bd prevalence that we observed.In autumn, less sunlight reaches canyon-bound streams, resulting in more limited basking opportunities [94,95], which may contribute to increased Bd prevalence in October.Basking opportunities are also dependent on microclimatic conditions, which are mediated by riparian vegetative cover and by flow regulation.Flood suppression by upstream dams allows encroachment of trees into stream channels, whereas wildfire can open the canopy.The Feather River watershed is characterized by a unique combination of dammed and free-flowing streams with frogs present in an area that has experienced intense wildfires [96].Rana boylii is designated as threatened under the US Endangered Species Act [36] in the Feather watershed, and frogs occur at low population densities relative to other regions [97].Ongoing population monitoring in this region will provide important insights into the relationships between Bd, flow regulation and canopy cover.
On the population level, persistence of R. boylii may be attributed to demographic resilience.Population growth rates of R. boylii positively correlate with summer stream temperature and are greater following wet years with high total annual volume of stream-flow [97].The longevity (greater than 10 years) and fecundity of some individuals [30,98], should allow populations to rebound from disease-related losses when instream conditions are favourable for recruitment.At the extreme, one Bd+ female swabbed in 2016 in the Pit River of northern California was recaptured in 2021 (K Breedveld 2021, personal observation).

Limitations and future work
We note some potential limitations of our study.Bd infections increased over time, but lower Bd prevalence in older museum samples may in part be attributed to lower detection probability.DNAdamaging formalin fixation compounded with degradation over time can limit DNA amplification, reducing qPCR detection probability relative to more recently preserved specimens and field-collected swabs [99].Thus, we cannot conclude definitively that prevalence was lower before 2005, or that the oldest Bd− museum samples truly indicate Bd's absence at a given location.Nonetheless, our samples royalsocietypublishing.org/journal/rsos R. Soc.Open Sci.11: 231270 were collected in similar time periods and regions as multi-species retrospective studies in the region that did detect Bd, and our sample sizes from earlier decades in areas where Bd was not detected in R. boylii are relatively robust.Thus, our detection data are likely to provide a reasonable estimate of Bd incidence across R. boylii populations over the last century.Another complication in estimating true prevalence is that Bd+ animals may have higher capture probabilities [69,100]; however this bias should be consistent between museum and field samples.In addition, our samples were not a random sample over the annual cycle of high and low flows in rivers that dictate frog activity: the majority of field samples in our dataset were collected during the dry season (May-October), which is the most biologically active period for R. boylii, when rivers can safely be waded and frogs can be encountered and caught.
We also note that SaTScan may not be optimal for river-dwelling organisms, as the linear distancebased method may cluster populations not connected by waterways.Sample sizes were uneven across regions and years, potentially causing some level of ascertainment bias; indeed, the contemporary clusters we detected are areas with fairly high R. boylii abundance and robust monitoring.However, a well-sampled and frog-abundant watershed (SF Eel River) was not part of a cluster due to a nonsignificant calculated relative risk.Finally, Bd infections do not necessarily indicate susceptibility to chytridiomycosis.For example, high pathogen loads can result from high vulnerability to Bd or high tolerance, i.e. the ability to maintain infection while remaining healthy.Assessing susceptibility requires coupling Bd assays with pathology, mortality/survivorship, or other disease indicators.
Our results suggest several avenues for future work to address Bd risk in R. boylii.Research that teases apart relationships between Bd and summer conditions (i.e.warm temperatures associated with reduced infections versus droughts associated with increased frog density and Bd transmission) will be especially informative for conservation.For example, stream restoration projects that aim to connect aquatic habitats could be designed to facilitate reduced disease transmission by considering frog density and thermal conditions, and reintroductions of zoo-reared frogs could be performed at a time of year that minimizes Bd infection risk.As cool temperatures predicted higher Bd infection probability, increased monitoring effort in autumn and spring may improve detection of chytridiomycosis outbreaks.Additional research is also needed to relate Bd-positivity to capture probability in this species, and to relate Bd loads to disease susceptibility among populations, life stages and sexes.Given the dual needs for data on susceptibility and minimizing impacts in remnant populations, we suggest expanding mark-recapture studies into clades designated as 'Endangered' (South Coast and South Sierra) and combining these with non-lethal tests of susceptibility such as immunocompetence assays [101], mucosal immunity assays [102], characterization of the microbiome [103][104][105] and its relationships with temperature across environmental contexts [106,107] and immunogenetics studies [74,75,108].
In addition, Bd genotype data in streams and rivers across California and Oregon are lacking.The global panzootic lineage (Bd-GPL) has been the primary genotype detected in California [60,[109][110][111][112][113].However, Bd-Brazil was detected on three amphibians in Alameda and Santa Clara counties (C Briggs 2019, personal communication; Nayar 2019, unpublished data), where we observed the highest Bd prevalence and loads in R. boylii.Virulence can significantly vary among Bd genotypes [114], and plastic gene expression in Bd-GPL can alter virulence [115].Indeed, Bd-GPL strains collected in northern California exhibit wide variation in virulence [13].Expanded knowledge of Bd genotype distribution may shed light on the regional differences we observed, and direct future evaluation of Bd's threat to R. boylii.
Addressing conservation challenges across large geographical scales requires extensive, multidisciplinary collaboration, especially for aquatic species in arid regions where there are competing demands for water and river flows are often over-allocated [116].Environmental flow frameworks designed to promote river health by allocating water for ecological purposes [117] would be enhanced by incorporating information on disease.In California, dam operators that generate hydroelectric power, irrigation districts that extract water for agriculture, and municipal utilities that impound rivers to supply drinking water to cities all play important roles in the conservation of riverine biota.There is a need for central reporting and collation of data collected by biologists across sectors including resource agencies, utility companies operating in dammed rivers, and academic scientists and community researchers working in undammed rivers to facilitate rapid knowledge sharing and coordination of responses to disease outbreaks.We present our compilation and analysis of diffuse data as a model for collaborative research on wide-ranging imperilled species.

Figure 1 .
Figure 1.Distribution of R. boylii samples assayed for Bd infection.Diamonds show museum samples (collected 1897-2005), circles show field samples (2005-2021).Filled symbols indicate Bd-positive samples.(a) Sampling locations across California and Oregon, USA.Symbols overlap in some localities; see inset barplots for sample sizes.Rana boylii clades are outlined and labelled, with California Endangered Species Act status abbreviated in parentheses: SSC = Species of Special Concern, TH = Threatened, EN = Endangered.(b) Spatio-temporal spread of Bd detections.Symbol size indicates sample size at the HUC-12 (sub-watershed) level.Generalized additive model (GAM) of latitude∼capture year + sample size in Bd-positive samples is shown with black curved line (R 2 = 0.133).Photo of R. boylii in Napa County, CA by Marina De León.

Figure 2 .
Figure 2. Batrachochytrium dendrobatidis infection prevalence ( per cent of individuals infected) in foothill yellow-legged frogs (Rana boylii) by (a) decade across museum and field samples, (b) month among field samples only (2005-2021).Museum samples are shown in gold and field samples in green.Error bars are Clopper-Pearson 95% confidence intervals, and dot size is proportional to total sample size.Blue shaded area plot in (b) (right y-axis) shows range of monthly precipitation and mean monthly precipitation (dotted line) across all sample sites.

Figure 3 .
Figure 3. SaTScan clusters of significant relative risk of Bd infection in Rana boylii.Samples were aggregated at sub-watershed (HUC-12) level.Historical clusters (a) are based on museum samples (collected 1897-2005); contemporary clusters (b) are based on field samples (2005-2021).Clusters are shown as large circles, with coloured dots at the centre of sub-watersheds included within a cluster, and small grey dots depicting underlying distribution of aggregated samples.High Bd risk (high-rate) clusters are shown in red (historical clusters 1-3, contemporary clusters 1, 3-6); low Bd risk (low-rate) clusters are shown in blue (contemporary clusters 2, 7).High-rate clusters are defined as areas with significantly higher relative risk than expected (RR > 1), while low-rate clusters have significantly lower relative risk than expected (RR < 1).

100 Figure 4 .
Figure 4. Influence of aridity and air temperature on Bd prevalence in R. boylii.Aridity index (SPEI) includes 24 months preceding frog capture; maximum air temperature (T max ) is the maximum over 30 days preceding capture.Asterisks indicate p < 0.05.(a-b) Marginal effects of aridity (a) and T max (b) on Bd detection probability predicted by GLMM.Red lines and dashed CIs indicate effects in high-rate clusters (cluster rate = 2; figure 3b); blue lines and CIs indicate effects in low-rate clusters (rate = 1; figure 3b); grey lines and CIs indicate effects in watersheds outside clusters (rate = 0).(c-f ) Linear regressions of aridity (c-d) or T max (e-f ) with observed Bd prevalence (% individuals infected) between rain-driven coastal systems versus inland snowmelt-driven Sierran systems.

Figure 5 .
Figure 5. Influence of aridity, body length and air temperature on Bd infection intensity in R. boylii.Aridity index (SPEI) includes 24 months preceding frog capture; frog length (SVL) was reported as snout-vent length; maximum air temperature (T max ) was the maximum over 30 days preceding capture.(a-c) Marginal effects of aridity (a), SVL (b), and T max (c) on Bd loads predicted by GLMM.Grey dashed lines show 95% CIs.Bd loads are reported as log-transformed ITS copy number.Linear regressions of aridity (d,e), SVL ( f,g), and T max (h,i) with observed Bd load between rain-driven coastal systems versus inland snowmeltdriven Sierran systems.